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ABSTRACT 

We use Bayesian model selection tools to forecast the Planck satellite's ability to distinguish 
between different models for the reionization history of the Universe, using the large an- 
gular scale signal in the cosmic microwave background polarization spectrum. We find that 
Planck is not expected to be able to distinguish between an instantaneous reionization model 
and a two-parameter smooth reionization model, except for extreme values of the additional 
reionization parameter. If it cannot, then it will be unable to distinguish between different 
two-parameter models either. However, Bayesian model averaging will be needed to obtain 
unbiased estimates of the optical depth to reionization. We also generalize our results to a 
hypothetical future cosmic variance limited microwave anisotropy survey, where the outlook 
is more optimistic. 
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1 INTRODUCTION 

The five-year data from the Wilkinson Microwave Anisotropy 
Probe (WMAP: Hinshaw et al. 2008; Dunkley et al. 2008; Ko- 
matsu et al. 2008) have given reasonably tight constraints on the 
optical depth to Thomson scattering from the last-scattering sur- 
face, r = 0.09 ± 0.02 with modest dependence on inclusion of 
additional datasets and changes to model assumptions. It has not 
however had the accuracy needed to go beyond this one-parameter 
description of the ionization history of the Universe, to give a more 
detailed view of how reionization took place and to distinguish be- 
tween the various models in the literature (though combined with 
tentative indication of a change in Lyman-a optical depth around 
redshift 7, it does give some indication that reionization is an ex- 
tended process). 

Theoretical studies suggest that the process of reionization can 
be quite complex (eg Barkana & Loeb 2001; Haiman & Holder 
2003; Cen 2003; for reviews see Barkana & Loeb 2007a; Meiksin 
2007). The Planck satellite may have the sensitivity to go beyond 
a one-parameter description of the process. For instance, Lewis, 
Weller & Battye (2006) considered three specific reionization histo- 
ries (with other cosmological parameters held fixed), and assessed 
whether Planck would be able to distinguish amongst them, finding 
that it did indeed have some ability to do so. 

However, the true data analysis problem is more complicated 
than in their study. Future experiments will not be trying to distin- 
guish between a small set of specific reionization histories. Rather, 
there will be competing models for reionization each of which fea- 
ture parameters that need to be determined from the data. That is, 
the problem is one of model selection (see Gregory 2005; Lid- 
dle, Mukherjee & Parkinson 2006a; Trotta 2008, and references 
therein). In this paper we use Bayesian model selection tools to 
forecast the ability of the Planck satellite, and a putative cosmic- 



variance-limited future survey, to distinguish between two reion- 
ization models, instantaneous reionization (parameterized solely by 
the optical depth to reionization r or equivalently the redshift of 
reionization), and a smooth transition to the ionized state, param- 
eterized by a further parameter d n which measures the rapidity of 
the transition (in conformal time rj). 

We consider only the large-scale bump in the cosmic mi- 
crowave background (CMB) polarization spectrum generated by 
Thomson scattering of the CMB quadrupolar anisotropy during 
reionization. The detailed shape of the bump is related to the evo- 
lution of the globally-averaged ionized fraction during reionization 
(Kaplinghat et al. 2003; Hu & Holder 2003; Colombo et al. 2005). 
The power on scales smaller than the horizon size at reionization 
is uniformly damped by e _2r ; this then cannot be used to con- 
strain the details of reionization beyond r, or even to constrain r 
itself which would be almost completely degenerate with the am- 
plitude of perturbations. Other degeneracies are discussed in Mar- 
tins et al. (2004) and Trotta & Hansen (2004). Reionization af- 
fects the CMB spectrum again on much smaller scales via sec- 
ondary effects due to inhomogeneous or patchy reionization and 
the Ostriker-Vishniac effects (Ostriker & Vishniac 1986; Weller 
1999; Hu 2000). We do not consider these effects which are be- 
yond multipole I ~ 2000, modelling only uniform reionization. 
In the future, 21cm emission from neutral hydrogen is expected to 
provide a good tracer of the details of reionization (eg. see Barkana 
& Loeb 2007b), and there are experiments that will focus on map- 
ping this emission. 



2 THE MODELS 

Our cosmological model is the usual spatially-flat ACDM cosmol- 
ogy, seeded by power-law adiabatic density perturbations. Its ad- 



© 0000 RAS 



2 Pia Mukherjee and Andrew R. Liddle 




justable parameters are the dark matter and baryon densities Q c 
and fib, the Hubble parameter h, and the perturbation amplitude 
A s and spectral index n s . These are fixed to WMAP3 best-fit val- 
ue£| (Spergel et al. 2007) for Q,\,h 2 , Q c h 2 , the projected sound 
horizon 9, A a exp(— 2r) and n s . We then study the reionization 
signal from the TE and EE spectra out to i of 100. It is possible 
to use such an analysis procedure because the non-reionization pa- 
rameters are very well determined by the TT spectrum, and because 
the large-scale signal in CMB polarization is independent of the 
other parameters. A similar procedure has been followed in works 
including Kaplinghat et al. (2003), Holder et al. (2003), and Mor- 
tonson & Hu (2008a,b). The uncertainty on r derived holding these 
parameters fixed is expected to be an underestimate by about 10% 
(Mortonson & Hu 2008a). 

We assume standard recombination. If the recombination 
model eventually needs to be modified to account for two-photon 
decays (Dubrovich & Grachev 2005; Wong & Scott 2007; Chluba 
& Sunyaev 2008; Hirata 2008), this should not affect the model 
comparisons we present here because it would be common to all 
the models. In addition, the spectrum changes on intermediate to 
small scales while we are using only the large scales here. 

We mainly consider a two-parameter reionization model de- 
fined by the ionization fraction history 



tanh 



XVi) = [1 -%e(Vi-l)] 



—2- - 1 



eL + l 



-HEe(f?i-l).(l) 



where x e refers to the ionization fraction, r\i and r\i-\ refer to con- 
secutive time steps, r\ to the conformal time at the i-th time step, 
z r is the redshift at which the ionization fraction is 0.5, rj Zl is the 
conformal time corresponding to that redshift, and d v gives the (in- 
verse) width of the transition. Such a transition is implemented in 
CAMB (Lewis, Challinor & Lasenby 2000), and the commonly- 
used instantaneous reionization scenario corresponds to d v having 
a large enough value, such as 50, that z T is effectively the redshift 
of instantaneous reionizationQ 

We additionally force the ionization fraction to unity for 
z < 6, to avoid conflict with quasar absorption spectrum data, and 
to zero for z > 30 as no ionizing sources are expected so early. 
The optical depth r is computed numerically for any such reioniza- 
tion history. Given the reionization history, we compute the CMB 
power spectra using a version of CAMB with minor modifications. 
Our assumed cosmological model has only scalar initial perturba- 
tions, and so we do not compute the BB polarization spectra. 

Figures Q] and [2] show some predicted power spectra for these 
models, showing in particular that r is indeed mainly responsible 
for variations in the predictions and hence the most readily mea- 
sured parameter. Furthermore, at fixed r it is clear from figure [2] 
that all the discriminating power is in the polarization spectra rather 
than temperature. 



2 We use a version of CAMB prior to April 2008, that includes only hy- 
drogen reionization and thus a final ionization fraction of unity. Including 



1 Our calculations predated the WMAP five-year data release (Hinshaw et 
al. 2008), which however left the numbers almost unchanged. 



helium reionization the final ionization fraction would be ~ 1.08. See sec- 
tion XIII of CAMB notes, via a link from http://camb.info/readme.html for 
further details. The model selection results presented in this paper are not 
expected to change following this inclusion of helium reionization. 
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Figure 2. As Figure[T] but for models each with r fixed at 0.1. At fixed r the spectra, especially TT, have only a weak dependence on 



Through most of the following analysis we take a fiducial 
z T -d v model with d v = 3 and z r = 8.9, corresponding to r = 0.1 
as already well determined by WMAP5. 



Flat priors are assumed on z r and d v , over ranges of 6-30 and 
0-10 respectively (the figures show that values of d v larger than 
this result in scenarios very close to instantaneous). Figure Q] also 
motivates us to try a logarithmic prior on d v , for which we take the 
range 0.3-30. However a r < 0.3 prior is also imposed, so that the 
one additional parameter as compared to instant reionization does 
not correspond to one additional degree of freedom. For this reason 
regular likelihood ratio tests, of the kind performed in Kaplinghat 
et al. (2003) and Holder et al. (2003), will not be valid. Figure [3] 
shows the induced prior on r resulting from linear priors on z T and 
d v (their own priors are not perfectly flat due to the extra imposition 
of the t < 0.3 prior). The uncertainty on r from Planck is of order 
0.01, and over such widths the prior on r is roughly uniform; it is 
even more so around the fiducial value of r. The same applies to 
the case of a log prior on d v . 



Other models of reionization have been proposed (Haiman & 
Holder 2003; Cen 2003), for instance the double reionization sce- 
nario considered in Lewis et al. (2006), which would in general re- 
quire a third parameter. From a model selection point of view (i.e. 
taking into account parameter uncertainties within each model) re- 
sults in this paper indicate that discriminating such a model from 
a smooth transition model would be beyond the scope of Planck, 
though perhaps within the scope of a closer to cosmic variance lim- 
ited experiment. Here our focus is mainly on clarifying what Planck 
can learn about reionization. 




Figure 3. Uniform priors on z T between 6 and 30, and on d r/ between and 
10, result in a non-uniform prior on r (solid curves). We can work with such 
a prior on r because expected uncertainties from a Planck-like experiment 
are At = 0.01 and over such a range the prior is fairly flat. 



3 MODEL SELECTION FORECASTING 
METHODOLOGY 

Model selection forecasting assesses a given experiment's ability 
to distinguish between different cosmological models. This ability 
necessarily depends on the true model and on its parameter values 
(and of course on the overriding assumption that one of the models 



© 0000 RAS, MNRAS 000. 000-000 



4 Pia Mukherjee and Andrew R. Liddle 





Model, priors 




parameter estimates 


In Evidence 


A In E 




instantaneous reionization 

z r : 6-30, dr, = 50 




Zr = 12.9 ±0.5 
t = 0.108 ±0.006 


-6.3 ±0.1 


0.0 


Planck satellite 


linear model 

z r : 6-30, d v : 0-10 




= 10.1 ±1.7, dr, =4.4 ±1.9 
t = 0.103 ±0.006 


-4.4 ±0.2 


1.9 




log d v model 
z r : 6-30, log dr,: -1.2-3.4 


Zr 


= 9.9 ±1.9, d^ =4.1 ±1.9 
r = 0.103 ±0.006 


-4.7 ±0.2 


1.6 




instantaneous reionization 

z r : 6-30, d v = 50 




z T = 12.7 ±0.2 
r = 0.106 ±0.003 


-15.8 ±0.1 


0.0 


cosmic variance 
limited 


linear dr, model 

z T : 6-30, dr,: 0-10 


Zr 


= 9.4 ± 1.0, d v = 3.3 ±0.5 
r = 0.102 ±0.005 


-6.6 ±0.1 


9.2 




log dr, model 
z r : 6-30, log dr,: -1.2-3.4 


Zr 


= 9.2 ± 1.0, d v = 3.2 ±0.4 
r = 0.101 ±0.005 


-6.8 ±0.3 


9.0 



Table 1. Analyzing TE and EE spectra of Planck specifications (first panel), with a fiducial model of z r =8.9 and dr, = 3 (implying r = 0.1) using three test 
models. The second panel shows the same for a cosmic variance limited experiment. In Evidences are based on four estimates of the evidence for each model. 



we will be considering is, if not the actual true model, at least repre- 
sentative of its predictions for the experiment under consideration). 
We call this true model and its parameter values the fiducial model. 
Trotta (2007) introduced an approach to model selection forecast- 
ing, PPOD, which averaged the model selection forecast over the 
current knowledge of model parameters, so as to give the prob- 
ability of the future experiment giving different model selection 
outcomes. Mukherjee et al. (2006b) adopted a different approach 
where the model selection outcome was forecasted as a function of 
the fiducial model parameters, so as to assess where in the parame- 
ter space the experiment could strongly distinguish the models, and 
defined some experimental figures of merit based on this notion. 

Computational restrictions prevent us from making an exten- 
sive investigation of how the model selection outcome will depend 
on the fiducial model chosen. Instead, we consider a single fidu- 
cial Zr-dr, model with d,, — 3 and z r = 8.9, corresponding to 
t = 0.1, and simulate Planck quality data for it. Figure [2] shows 
that this model is quite different from an instantaneous reioniza- 
tion model with the same r (though models can be constructed that 
are even more drastically different, up to the step model seen in 
Figure 1). Our goal is to determine whether Planck could distin- 
guish the two-parameter model from an instantaneous reionization 
model, for these fiducial parameters with a certain level of strength 
of evidence. If it can, then the threshold for detection lies between 
this model and the instantaneous reionization model, otherwise it 
lies further away. 

Throughout we use the nested sampling algorithm for comput- 
ing model evidences, first implemented for cosmological applica- 
tions in Mukherjee, Parkinson & Liddle (2006a) and subsequently 
developed in Parkinson, Mukherjee & Liddle (2006). It was used to 
make forecasts for Planck's ability to determine inflationary param- 
eters in Pahud et al. (2006,2007). This algorithm, due to Skilling 
(2006), computes the Bayesian evidence for any given model, as 
well as providing parameter estimates within that model. The evi- 
dence is the probability of the data given the model, hence can be 
used to determine how likely each model is to have given rise to the 
data. A difference of 2.5 in log evidence can be taken to be signif- 
icant, and 5 decisive, evidence in favour of the model with larger 
evidence (Jeffreys 1961). 



4 RESULTS 

4.1 The Planck satellite 

We model Planck TE and EE data using just the 143 GHz polar- 
ization channel, following for its specifications the current Planck 
documentation^ The full likelihood is constructed in the manner of 
Lewis (2005) and Pahud et al. (2006,2007), without creating noisy 
data realizatons. This ensures that the bias issues we discuss be- 
low are not a result of realization noise, but instead the forecast is 
equivalent to averaging over many data realizations and is thus it- 
self effectively unbiased (see the appendix of Sahlen et al. 2008). 
We assume a sky coverage of 0.8, and take the likelihood up to a 
maximum multipole of £ = 100. 

The first entry in Table [T] shows that the result of analyz- 
ing Planck data based on this chosen fiducial model with a one- 
parameter instantaneous reionization model (ie. with d n held fixed 
at 50, varying z T over the 6-30 prior range). The second and third 
entries show the data analyzed with the (correct) Zr-d v model. Lin- 
ear and log priors on d v are assumed, to check for the dependence 
of results on such assumptions. 

Our main result is the relative evidences of these models, 
where the instantaneous reionization model has a In evidence which 
is less than 2 smaller than the two-parameter models. Accordingly, 
Planck will not be good enough to exclude the instantaneous reion- 
ization model, even though the true model appears to have a quite 
different ionization history. Put another way, Planck is not powerful 
enough to explore two-parameter models of reionization (at least 
unless the true model is even further from instantaneous reioniza- 
tion than our fiducial model). 

Besides the evidence, one can also compute the Bayesian com- 
plexity of Planck data for the chosen fiducial model in the manner 
of Kunz, Trotta & Parkinson (2006). Using such an analysis Kunz 
et al. found that r is already a required parameter with WMAP 3-yr 
data (and a similar analysis would likely show that another reion- 
ization parameter is required when the evidence supports it). 

This is, however, not quite the end of the story. The parameter 
distributions given from the nested sampling algorithm are shown 
in Figure[4] It is apparent from this that the estimated r (and z T ) are 

3 www.rssd.esa.int/index.php?project=PLANCK&page=perf_top 
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Figure 4. A fiducial model with z T = 8.9 and = 3 treated with an instantaneous reionization model fixed at 50, z r allowed to vary between 6 and 30) 
(left panel), a z r -d v reionization model with z r between 6 and 30 and d v between and 10 (centre panel), and a z r -d v reionization model with log prior in 
dr, (z r between 6-30, log d v between -1.2 and 3.4 corresponding to 0.3 and 30) (right panel). 




0.095 0.1 0.105 0.11 0.115 0.12 0.125 
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Figure 5. As Figure|4] for a cosmic variance limited experiment. The axes scales are the same as in Figure|4] 



biased high in the instantaneous reionization model. The bias goes 
away in the two-parameter model, as it should since that model 
can describe the true behaviour of the data. Accordingly, to avoid 
a possible bias in measuring r one should consider both the one- 
parameter and two-parameter models, and Bayesian model average 
as in Liddle et al. (2006b) to obtain constraints on r (the cost be- 
ing a slightly increased uncertainty in r). Similar conclusions have 
been reported in other papers (eg. Kaplinghat et al. 2003; Holder et 
al. 2003). 

We have made some assumptions about the true (fiducial) 
model that we are not certain about. In practice we don't know 
the fiducial jz max when reionization started. This has been assumed 
to be 30 in the fiducial model. A different z m ax corresponds to a 
different reionization history, hence z ma x could be treated as an ad- 
ditional reionization parameter, but we don't go into a third reion- 
ization parameter here. Instead we ask what outcome arises if we 
analyze data so simulated (with a z ma x of 30) with a z T -d v model 
with z m ax = 20. Such an 'incorrect' model would not be distin- 
guishable from the true model by Planck. Further, if our incorrect 
model was not a smooth transition model but one involving a step 



function, again with only two parameters, corresponding to z max 
(prior range 7-30), with reionization ending at redshift 6, and with 
a constant reionization fraction in between these two redshifts of x e 
(prior range 0-1), such a model would again not be distinguishable 
from the assumed true model by Planck. These results are borne out 
of numbers presented in the next subsection for a cosmic variance 
limited hypothetical experiment. 



4.2 Cosmic variance limited case 

For a cosmic variance limited hypothetical experiment, the corre- 
sponding results are shown in the lower panel of Table Q] and in 
Figure [5] Again only TE and EE spectra are considered out to a 
maximum multipole of 100. This time the evidence favours the 
smooth and gradual transition z T -d v model decisively over the in- 
stantaneous reionization model. 

These results also show that, as before, a simpler model leads 
to a biased r, a bias that disappears upon using a complicated 
enough model for reionization. The choice of prior on d v (log or 
linear) doesn't make much difference. 
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Model, priors parameter estimates In Evidence A In E 



8.8 ± 1.4, d v = 1.7 ± 0.6 -8.0 ±0.1 7.8 

t = 0.100 ±0.004 

18.9 ± 0.8, x e = 0.39 ± 0.05 -10.8 ±0.1 5.0 
t = 0.094 ± 0.003 



Table 2. As Table [T] but analyzing the same fiducial model with an incorrect test model, for the cosmic variance limited case. The A In E are with respect to 
the instantaneous model in Table[T] 



(incorrect) z max = 20 z r = 

z r : 6-20, d v : 0-10 
(incorrect) z max -£ e model z max = 
Zmax: 7-30, x e : 0-1 



Table [2] shows an incorrect assumption regarding z max of 
the Zr-dq model does not significantly affect the evidence, ie. 
Zmax = 20 is indistinguishable from z m ax = 30. d v is underes- 
timated to make up for the difference in 2 max , and r is not misesti- 
mated. It also shows that the smooth and gradual transition model is 
favoured strongly, almost decisively, over the incorrect step model 
based on the difference in log evidence, while r is biased low under 
this incorrect model assumption. Both incorrect models are clearly 
distinguishable from (and favoured over) the instantaneous reion- 
ization model. However the fact that different choices of z max are 
not distinguishable indicates that even cosmic variance limited ex- 
periments cannot probe very fine details of the reionization history. 



5 CONCLUSIONS 

We find that Planck is not expected to be able to distinguish signif- 
icantly between a single-parameter and a two-parameter model of 
reionization in the model comparison sense, though it will mildly 
favour the two-parameter model for our chosen fiducial values. 
If the parameter values of the two-parameter true model were 
more extreme, then Planck might favour it significantly. However 
Bayesian model averaging the parameters over the two models will 
eliminate any bias in the optical depth to rescattering. 

A cosmic variance limited hypothetical experiment will be 
able to decisively distinguish between the one- and two-parameter 
models, to distinguish between some two-parameter models, and 
may be able to go onto a third parameter. 

The model comparison approach advocated here should if 
possible be applied to models parameterized by describing physi- 
cally relevant quantities for the reionization history model; the phe- 
nomenological quantities employed here are an intermediate step 
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